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ABSTRACT 

Excessive costs prohibit the Department of Defense and 
industry from testing entire systems to obtain estimates of 
system reliability. Thus, individual component tests must 
be combined in such a way as to yield an accurate estimate 
of the system reliability. 

Assuming each component reliability to have a prior 
distribution which is beta distributed, Bayesian techniques 
result in the posterior distributions also being beta dis- 
tributed. The system reliability for a series system would 
then be the product of the posterior distributions. 

Computer simulation techniques were used to determine 
the system reliability distribution for a series system of 
beta distributed components. Method of Moments techniques 
were then used to fit beta and normal distributions to the 
system distribution. The twentieth percentile points of 
the fitted and the actual distributions were then compared 
as a measure of accuracy of the fit. 

The fit of the beta distribution proved to be accurate 
for all parameter ranges and number of components in each 
system. The fit of the normal distribution was accurate 
only when used with limited parameter values. 
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I. INTRODUCTION AND SUMMARY 



A. INTRODUCTION 

The purpose of this study was to analyze the accuracy 
of fitting a product of beta distributions with beta and 
with normal distributions using computer simulation tech- 
niques. If successful, a normal or beta distribution could 
be used to fit the posterior distributions in Bayesian 
techniques to the reliability of a series system where the 
component prior distributions are also beta distributed. 

B. MEASURE OF ACCURACY 

The Bayesian 100(l-a)% lower confidence limit for the 
system reliability, R , is the 100 ath percentile point of 
the distribution of R . Therefore, the difference between 
the 100 ath percentile points of the actual and the fitted 
distributions was selected as a measure of the accuracy of 
the fit. 

C. SUMMARY OF RESULTS 

When the method of moments technique was used to fit a 
normal distribution to the distribution of the reliability 
of a series system, the fit was accurate for a limited num- 
ber of cases. The data depicted in Table III indicates 
that the normal distribution fit can offer a simple and ac- 
curate method of estimating the system reliability. 

When the constraints of the fitted normal cannot be 
met or if an easy method of obtaining a percentile point 
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of a beta distribution is available, the beta fit should be 
used. 
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II. THE PROBLEM 



A. BACKGROUND 

Due to the costliness of mission testing entire sys- 
tems, the Department of Defense must rely on individual 
component tests to estimate system reliability. Among 
other things, the accuracy of this estimate is dependent 
upon the fit of a specific distribution to actual distri- 
butions. It is often necessary due to mathematical and/or 
statistical complications to assume a distribution which 
only approximates the actual distribution. A secondary 
yet important consideration is the manner in which the com- 
ponent test data is utilized in determining the estimated 
system reliability. 

Bayesian techniques offer an acceptable solution to 
the test data utilization. In Bayesian analysis, the beta 
distribution with parameters a and b is usually choosen as 
the prior distribution for component reliability. The pa- 
rameters a and b relate to attribute type data. A problem 
arises when trying to combine the component beta distribu- 
tions to get a confidence interval on system reliability. 

In particular, the problem exists on how to construct methods 
for estimating the system reliability when each of the com- 
ponent reliabilities is assumed to be beta distributed. 

B. BETA DISTRIBUTION AND BAYESIAN TECHNIQUES 

1 . Beta Distributions 

The beta distribution with parameters a and b has 

the following probability function: 
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a ,b ,x>0 . 



ffy-) = T (a+b) x a " 1 1 

UXJ r(a)T(b) X U XJ 

The mean, y, and variance, o 2 , are computed by the following 
formulae : 

a 

y = I^b 

a 2 = ab 

(a+b) (a+b+1) 

2 . Bayesian Techniques 

The usefulness of beta prior distributions and 
Bayesian techniques can be depicted by the following example 
Suppose n components of one type are mission tested and 
further suppose each component is beta distributed, 8(x;a,b) 
Let s be the number of components which do not fail. The 
100(l-a)% lower Bayesian confidence limit on component reli- 
ability is the solution for P in the equation 

P 

a = / $(x;a+s ,a+b+n-s+l)dx. 

0 

That is, the posterior distribution of component reliability 
is 3 (x ;a+s ,b+n-s) . Thus, by the use of beta distributions 
and Bayesian techniques, both past and present test data can 
readily utilized in reliability calculations. 
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III. SOLUTION 



A. ACTUAL SYSTEM DISTRIBUTION 

The actual system distribution to which the beta and 
normal distributions were fitted was generated on the IBM 
360 Computer at the United States Naval Postgraduate School, 
Monterey, California. The main program was used with slight 
modifications to evaluate the lOOath percentile points of 
the fitted beta distributions. The- required subroutines 
were obtained from the Computer Center Library. 

1 . Main Program 

The main program was constructed to randomly gen- 
erate a reliability value for each component of the series 
systems. Assuming the component distributions to be inde- 
pendent, the system reliability, R , would then be the pro- 

S K 

duct of the component reliabilities, R^ ; i.e., R s =^n^ R^ . 

Thus, by multiplying the generated component reliability 
values, a value of the system reliability was obtained. 

This process was repeated to obtain 500 values of R g . The 
values were ordered to obtain a simulated distribution of 
R g . The one hundredeth number of this distribution was then 
assumed to be the twentieth percentile point of the actual 
system reliability distribution. 

2 . Subroutines 

The Subroutine BDTR was used to obtain values of 
the cumulative distribution function of the beta distributions. 
The subroutine is contained in R. E. Bargmann and S. P. Ghosh, 
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Statistical Distribution Programs for a Computer Language, 
I.B.M. Research Report RC-1094. It has a maximum error of 

- 5 

10 for values of a and b between 0.5 and 10 . Subroutine 
BDTR used subroutine NDTR, CDTR, and DLGAM in its .computa- 
tions . 



B. METHOD OF MOMENTS BETA FIT 

The supposition is that the system reliability, R , is 

K S 

beta distributed with parameters A and B where R =.n, R. 

r s 1=1 1 

and the R.'s are beta distributed with parameters a- and b. 
respectively. Using the method of moments technique, the 
mean, y, and variance, a 2 , of the system distribution are 
as follows: 



y = E(R s ) = J 1 E(R.) 



A 

A+B 



K a i ( a i + l) 

i = l (a.+b.Ha.+b.+l) ' 



o 



2 



AB 

(A+B) z (A+B+l) 



E(R S ) 2 



[E(R S )] 2 . 



In addition, 



Letting 



and 



2 , ACA+l) k a.Ca^l) 

E ^ R s ^ " (A+B) (A+B+l) " iSl (a i +b i ) (a i +b i +l) ’ 



K a i (a i +l) 

^ 2 ~ i^l (a- +b . ) (a- +b . +1) ’ 

^ i i J v l l 
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A and B can be determined: 



a = Cdi-ya) u_i 
y 2 -y i 2 



r = Ci-yi) (yi -y2) 
y 2 -yi 2 

Subroutine BDTR was then used to determine the twentieth 
percentile point of the beta distribution with parameters A 
and B as calculated. 



C. METHOD OF MOMENTS NORMAL FIT 

Define Q.=l-R. and assume R. is large." Q. is a beta 

x i l l 6 x i 

distributed random variable with parameters b. and a.. R 
K i i s 

is equal to ^n^(l-Q^). Since is "small," expanding 

K K 

i g 1 (l-Q i ) and neglecting the cross products, R s = l-^I|-LQi and 

K 

l'Rs = i2l Q i* The su PP os i t i° n as now that 1-R is normally 
distributed, with parameter y and a 2 since it is a sum of 

random variables. Using method of moments techniques, 

K K Kb- 

y = E (1-R )=E(.Jl, Q.) = .n. E(Q.) = .n 1 — -±- 
v s J v i=l x i 2 i=l ^ x i' i=l a^+b^ 



and 



K K 

a 2 = var(l-R s )=var( i 0 1 Qi) = i2 1 varCQ^ 



a . b • 

= n - 11 

i-1 (a- +b . ) 2 (a . +b . +1) 
v l i J K i i J 



The twentieth percentile point of R g is the same as the 
eightieth percentile point of 1-R s . The eightieth percentile 
point of 1-R^. is a Z Q +y where Z 0 is the eightieth percentile 
point of the Standard Normal Distribution. 
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IV. RESULTS 



A. ACTUAL SYSTEM DISTRIBUTION 

The actual system distributions were generated for 
series systems of fifteen and twenty-five components using 
the nineteen sets of parameters as listed in the first 
three tables. It was assumed that all components in the 
systems would have the same parameter values. The twentieth 
percentile points of these distributions are listed in 
Tables I and I I . 

Due to the number of runs required and to the length 
of computer time required to call subroutine BDTR, the num- 
ber of samples forming these distributions was limited to 
five hundred. This limitation resulted in the possible 
variance of the twentieth percentile point of a given dis- 
tribution by as much as .003 when a different series of ran- 
dom numbers were used to generate the same distribution. 
Variations of this amount would make a significant change 
in column five of Tables I and II. For example, using 
a=1000.0 and b=5.0 in Table I, the difference between the 
actual twentieth percentile point and the fitted beta point 
was only .002. Yet, this difference resulted in an error 
of 4.6 percentage points (20.0-15.4=4.6). 

B. METHOD OF MOMENTS BETA FIT 

Using the formulae developed previously to obtain the 
system parameters A and B, the twentieth percentile point 
of the system distribution was calculated by computer using 
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subroutine BDTR. Considering the maximum error of BDTR, 
the data of column four in Tables I and II should be more 
accurate than that of the actual distribution. 

Examination of columns three and four of Tabl.es I and 
II shows that the beta fitted distribution provides a more 
accurate estimate of the actual distribution for the system 
of twenty-five components than for the fifteen. However, 
even for the smaller system, the method does result in rea- 
sonably accurate estimates of the twentieth percentile point 
particularly if a variance of .003 is introduced into the 
actual distribution percentile point. It is interesting to 
note that the method appears to be about the same degree ac- 
curate throughout the parameter ranges for both systems. 

Due to the steepness of the curves involved, the small per- 
centage error figures can be misleading. The actual per- 
centile points (column five) might offer a better measure 
of accuracy. These points indicate that in most of the 
cases the error is conservative and that the approximation 
actually represents the lower limit of a higher degree of 
confidence than eighty percent. 

C. METHOD OF MOMENTS NORMAL FIT 

The results of the normal fit are depicted in Table 
III. As can readily be seen in this table, the normal fit 
is not accurate when the ratio of a to b is less than 500 
to 1. Since these parameters can be related to the number 
of failures and number of successes, these results tend to 
prove the initial assumption made when deriving this method. 
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The results in- 



This assumption was that was "large." 
dicate that "large" means one or less failures out of 500 
mission tested components. Given these parameter con- 
straints, the normal fit is quite accurate. Due to the 
limited range of parameters, this method was not investi- 
gated further. 

D. BETA FIT WITH VARIOUS PARAMETERS 

The process of fitting the beta distribution to the 
system reliability was successful enough to investigate the 
accuracy of the method when the components have different 
parameter values. Table IV depicts these results. The 
comparison was made for systems of fifteen and twenty-five 
components and for two different general parameter ranges. 

The limitations of computer core suggested that no more 
than five different distributions be used. Thus, in the 
fifteen component system, there are three components with 
each set of parameters; and in the twenty- five component sys- 
tem, there are five components with each set of parameters. 

The results indicate that the method of moment fitted 
beta distribution is equally accurate when there are mixed 
distributions in the systems. 
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TABLE I 



This table lists the results of the method of moments 
beta fit and the actual distribution for a series .system of 
fifteen components. All distributions in the system have 
the same parameters. Columns one and two list the a and b 
parameters respectively. Column three is the twentieth 
percentile point of the actual distribution. Column four 
is the twentieth percentile point of the beta fitted dis- 
tribution. Column five lists the percentile point of the 
actual distribution that the number in column four represents. 
Column six is the percent error which is the difference be- 



tween 


column 


three and column 


four divided 


by column 


three 


a 


b 


Actual 


Beta Fit 


Act . % 


Error 


10.0 


0.5 


0.372 


0.372 


20.0 


0.0 


10.0 


1.0 


0.168 


0.163 


18.8 


3.0 


30.0 


0.5 


0.768 


0.723 


4.4 


5.9 


30.0 


1.0 


0.545 


0.547 


20.8 


0.3 


30.0 


2.0 


0.328 


0.323 


17.4 


1.6 


50.0 


0.5 


0.825 


0.824 


19.4 


0.1 


50.0 


1.0 


0.700 


0.696 


18.8 


0.6 


50.0 


2.0 


0 . 529 


0.505 


7.8 


4.5 


100.0 


0.5 


0.913 


0.908 


16.4 


0.6 


100.0 


1.0 


0.839 


0.834 


17.6 


0.6 


100.0 


2.0 


0 . 716 


0.709 


16.4 


1.0 


500.0 


0.5 


0.981 


0.981 


20.1 


0.0 


500.0 


1.0 


0.965 


0.964 


17.8 


0.0 


500.0 


2.0 


0.936 


0.933 


13.0 


0.5 


500.0 


5.0 


0 . 850 


0.849 


17.4 


0.2 


1000.0 


0 . 5 


0.990 


0.990 


20.1 


0.0 


1000.0 


1.0 


0.983 


0.982 


17.6 


0.0 


1000.0 


2.0 


0.967 


0.966 


17.6 


0.0 


1000.0 


5.0 


0.923 


0.921 


15.4 


0.2 
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TABLE II 



This table lists the results of the method of moments 
beta fit and the actual distribution for a series .system 
of twenty-five components. All distributions in the sys- 
tem have the same parameters. Columns one and two list the 
a and b parameters respectively. Column three is the twen- 
tieth percentile point of the actual distribution. Column 
four is the twentieth percentile point of the beta fitted 
distribution. Column five lists the percentile point of 
the actual distribution that the number in column four 
represents. Column six is the percent error which is the 
difference between columns three and four divided by column 
three . 



a 


b 


Actual 


Beta Fit 


Act . % 


Error 


10.0 


0.5 


0 . 205 


0 . 207 


20.8 


0.8 


10.0 


1.0 


0.055 


0.054 


19.4 


1.8 


30.0 


0.5 


0.650 


0 . 597 


2.6 


8.2 


30.0 


1.0 


0.379 


0.379 


20.0 


0.0 


30.0 


2.0 


0.160 


0.160 


20.1 


0.0 


50.0 


0 . 5 


0.733 


0.735 


20.8 


0.3 


50,0 


1.0 


0.560 


0.559 


19.6 


0.2 


50.0 


2.0 


0.358 


0.331 


5.0 


7.5 


100.0 


0 . 5 


0.860 


0.857 


17.6 


0.3 


100.0 


1.0 


0 . 747 


0.747 


20.0 


0.0 


100.0 


2.0 


0.575 


0.574 


19.4 


0 . 2 


500.0 


0.5 


0.970 


0.970 


20.8 


0.0 


500.0 


1.0 


0.944 


0.944 


19.6 


0.0 


500.0 


2.0 


0 . 897 


0. 894 


14.2 


0 . 3 


500.0 


5.0 


0.765 


0.765 


21.0 


0.0 


1000.0 


0 . 5 


0.985 


0.985 


21.0 


0.0 


1000.0 


1.0 


0.971 


0.971 


19.0 


0.0 


1000.0 


2.0 


0.945 


0.946 


20 . 8 


0.1 


1000.0 


5.0 


0.875 


0.875 


18.0 


0.0 
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TABLE III 



This table lists the results of the method of moments 
normal fit and the actual distribution for a series system 
of twenty- five components. All distributions in the system 
have the same parameters. Columns one and two list the a 
and b parameters respectively. Column three is the twen- 
tieth percentile point of the actual distribution. Column 
four is the twentieth percentile point using the techniques 
described previously. 



a 


b 


Actual 


Normal Fit 


10.0 


0.5 


0.205 


- 0.454 


10.0 


1.0 


0.055 


- 1.62 


30.0 


0.5 


0.650 


0.49 


30.0 


1.0 


0.379 


0.062 


30.0 


2.0 


0.160 


- 0.739 


50.0 


0.5 


0.733 


0.694 


50.0 


1.0 


0.560 


0.429 


50.0 


2.0 


0 . 358 


- 0.072 


100.0 


0 . 5 


0.860 


0.846 


100.0 


1.0 


0.747 


0.711 


100.0 


2.0 


0.575 


0.452 


500.0 


0.5 


0.970 


0.969 


500.0 


1.0 


0.944 


0.942 


500.0 


2.0 


0.897 


0.888 


500.0 


5.0 


0.765 


0.734 


1000.0 


0 . 5 


0.985 


0.984 


1000.0 


1.0 


0.971 


0.971 


1000.0 


2.0 


0.945 


0.944 


1000.0 


5.0 


0.875 


0.866 
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TABLE IV 



This table lists the results for method of moments beta 
fit using different parameters in the component distributions. 
Columns one and two are a and b parameters used. Column 
three is the number of components in the system. Five dif- 
ferent sets of parameters are used. In the fifteen compo- 
nent system, there are three components with the same set of 
parameters. In the twenty-five component system, there are 
five components with the same set of parameters. Column 
four is the twentieth percentile point of the actual distri- 
bution. Column five is the twentieth percentile point of 
the fitted beta distribution. 



a 


b 


k 


Actual 


Beta Fit 


50.0 


0.5 








100.0 


0.5 


15 


0.822 


0.820 


100.0 


1.0 








100.0 


2.0 


25 


0.739 


0.728 


500.0 


5.0 









500.0 

500.0 


1.0 

2.0 


15 


0 .966 


0.966 


1000.0 

1000.0 


0.5 

1.0 


25 


0.946 


0.947 


1000.0 


2.0 
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TABLE V 



Using the method of moments beta fit, R g was assumed to 
be beta distributed with parameters A and B. This, table 
lists the parameters A and B for the various numbers of sys- 
tem components and component parameters a and b. Columns 



one and 


two are the 


a and b parameters respectively 


. Col- 


umns three and six 


are the number of 


components in 


the sys 


tern. Columns 


four 


and five are the A 


and 


B values 


for the 


variables of 


columns one, two 


, and three. 


Columns 


seven 


and eight are 


the A 


and B values for 


the 


variab les 


of col- 


umns one 


! , two 


, and 


six . 










1 . a 


2 . b 


3 . k 


4. A 


5 . B 


6.k 


7. A 


8 . B 


10.0 


0.5 


15 


7.2 


7.5 


25 


5.9 


14.0 


10.0 


1.0 


15 


5.5 


17.5 


25 


3.8 


37.8 


30.0 


0.5 


15 


26.8 


7.5 


25 


24.8 


12.7 


30.0 


1.0 


15 


24.1 


15.3 


25 


20.8 


26.4 


30.0 


2.0 


15 


19.8 


32.3 


25 


15.3 


61.4 


50.0 


0.5 


15 


46.7 


7.5 


25 


44.5 


12.6 


50.0 


1.0 


15 


43.7 


15.1 


25 


39 . 8 


25.5 


50.0 


2.0 


15 


38.5 


30.9 


25 


32.4 


54.0 


100.0 


0.5 


15 


96.6 


7.5 


25 


94.3 


12.5 


100.0 


1.0 


15 


93.3 


15.0 


25 


88.9 


25.1 


100.0 


2.0 


15 


87.4 


30.2 


25 


79.6 


51.0 


500.0 


0.5 


15 


500 . 7 


7.6 


25 


498.9 


12.6 


500.0 


1.0 


15 


490.7 


14.9 


25 


485.4 


24.9 


500.0 


2.0 


15 


485.9 


30.0 


25 


476.1 


50.0 


500.0 


5.0 


15 


465.3 


74.9 


25 


443.4 


125.2 


1000.0 


0.5 


15 


1029.9 


7.8 


25 


1046.6 


13.2 


1000.0 


1.0 


15 


976.0 


14.7 


25 


969.2 


24.5 


1000.0 


2.0 


15 


993.8 


30.2 


25 


981.4 


50.3 


1000.0 


5.0 


15 


966.8 


75.1 


25 


942.9 


125.2 


50.0 


0.5 


3 






5 






100.0 


0.5 


3 






5 






100.0 


1.0 


3 


89.98 


16.03 


5 


85.27 


26.79 


100.0 


2.0 


3 






5 






500.0 


5.0 


3 






5 






500.0 


1.0 


3 






5 






500.0 


2.0 


3 






5 






1000.0 


0.5 


3 


578.16 


16.69 


5 


576.65 


28.02 


1000.0 


1.0 


3 






5 






1000.0 


2.0 


3 






5 
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MAIN PROGRAM FOR FIVE DIFFERENT COMPONENT TYPES 



DIMENSION RSHAT ( 500) ,XAR AY ( 5,750) , P ARAY ( 5 , 75 0 ) , NUM ( 5 ) 
Z=URN ( 0 ) 

DATA NITER, 137, I 3 , KCOMP , PROD , NRE P , ALPER , I 3 1 / 

DO 555 132=1,5 

400 READ( 5,401)ALFA,BETA,NSTAR,NST0P,N,CFACT 

401 FORMAT (2F6.1»2If» 1 1 , F 8 • 7 ) 

I 2=1 

DO 402 I21=NSTAR,NST0P,N 
X=I21*CFACT 

CALL BDTR( X, ALFA, BETA, P,DOUT , IER) 

XARAY ( 132, I 2 ) = X 
PARAYI 132, I2)=P 
12=12+1 

402 CONTINUE 
NUM( I 32 ) = 1 2-1 

555 CONTINUE 

DO 222 11=1, NITER 
132=1 

DO 333 14=1, KCOMP 
FLY=URN ( 1 ) 

49 IF(FYL-PARAY( 132,13) )50, 50, 51 
51 13=13+1 

I F ( I3.EQ.NUMU32 ) )G0 TO 60 
GO TO 49 

50 A=PARAY( 132,13 J-FYL 
I F( I 3« EQ« 1 ) GO TO 60 
I T3=I3- 1 

B=FYL~PAR AY ( 132, IT3 ) 

I F ( A— B )60, 60,61 

61 RIHAT=XARAY( 132, IT3 ) 

GO TO 62 

60 RIHAT=XARAY( 132, 13) 

62 PROD=PROD*RIHAT 
13=1 

A 14=14 
AM=AI4/I37 
MM= AM+ . 2 
M=I4/ 137 

IF(MM.EQ.M) 132=132+1 
333 CONTINUE 

RSHAT ( I 1 ) = PROD 
PROD= 1 * 

222 CONTINUE 
J=0 

8 IF( J. EQ. NITER-1 )G0 TO 11 
1 = 1 

J = 0 

9 IF (RSHATt I ) .GT.RSHATI 1+1) ) GO TO 10 
J=J+1 

1 = 1 + 1 

I F ( I .EQ.NITER)GO TO 8 
GO TO 9 

10 RTEMP=RSHAT( I) 

RSHAT ( I ) =RSHAT ( I +1 ) 

RSHAT ( 1+1 )=RTEMP 

1 = 1 + 1 

I F C I . EQ «NI TER ) GO TO 8 
GO TO 9 

11 CONTINUE 
STOP 

END 
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